##########################
### This file generates all the models used in the paper
##########################

###############################
### Dyadic-year main models ###
###############################

## DY LDV
coop.f.d.l <- log(post.coopTotal + logdelta)  ~  log(coopTotal + logdelta) + log(abs(hostTotal) + logdelta) +
  (KOFEcGIdf.x.0 + KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dep1 * dep2 + dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) +
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2+ 
  alliance + I(contiguity>0)

conf.f.d.l <- log(abs(post.hostTotal) + logdelta)  ~ log(coopTotal + logdelta) + log(abs(hostTotal) + logdelta) + 
  (KOFEcGIdf.x.0 + KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dep1 * dep2 + dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) +
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2+ 
  alliance + I(contiguity>0)

fitsur.f.d.l <- systemfit(list(coop = coop.f.d.l, conf = conf.f.d.l), data=final.data.dyad, method = "OLS")
summary(fitsur.f.d.l)

sur.f.d.l <- systemfit(list(coop = coop.f.d.l, conf = conf.f.d.l), data=final.data.dyad, method = "SUR")
summary(sur.f.d.l)

## Testing whether SUR is necessary
lrtest(fitsur.f.d.l,sur.f.d.l)


### present relevant linear combinations
linearHypothesis(fitsur.f.d.l, "coop_KOFEcGIdf.x.0 = conf_KOFEcGIdf.x.0", test = "Chisq")
linearHypothesis(fitsur.f.d.l, "coop_KOFEcGIdf.y.0 = conf_KOFEcGIdf.y.0", test = "Chisq")
linearHypothesis(fitsur.f.d.l, "coop_major.sym = conf_major.sym", test = "Chisq")
linearHypothesis(fitsur.f.d.l, "coop_major.sym + 0.117 * coop_KOFEcGIdf.x.0:major.sym = conf_major.sym + 0.117 * conf_KOFEcGIdf.x.0:major.sym", test = "Chisq")
linearHypothesis(fitsur.f.d.l, "coop_major.sym -0.128 * coop_KOFEcGIdf.x.0:major.sym = conf_major.sym -0.128 * conf_KOFEcGIdf.x.0:major.sym", test = "Chisq")


######################### -----------------------------------------------------
## Models in the appendix
#########################------------------------------------------------------

## 1. DY Year FE -----------------------------------------------------------------
coop.f.d <- log(post.coopTotal + logdelta)  ~   (KOFEcGIdf.x.0 + KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dep1 * dep2 + dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) + 
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2 + as.factor(year)+ 
  alliance + I(contiguity>0)

conf.f.d <- log(abs(post.hostTotal) + logdelta)  ~  (KOFEcGIdf.x.0 + KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dep1 * dep2 + dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) + 
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2 + as.factor(year)+ 
  alliance + I(contiguity>0)

fitsur.f.d <- systemfit(list(coop = coop.f.d, conf = conf.f.d), data=final.data.dyad, method = "OLS")
summary(fitsur.f.d)
## ----------------------------------------------------------------------------

## 2. DY LDV three way interaction --------------------------------------------
coop.f.d.l.i <- log(post.coopTotal + logdelta)  ~  log(coopTotal + logdelta) + log(abs(hostTotal) + logdelta) +
  (KOFEcGIdf.x.0*KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dep1 * dep2 + dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) +
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2 + 
  alliance + I(contiguity>0)

conf.f.d.l.i <- log(abs(post.hostTotal) + logdelta)  ~ log(coopTotal + logdelta) + log(abs(hostTotal) + logdelta) + 
  (KOFEcGIdf.x.0*KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dep1 * dep2 + dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) +
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2 + 
  alliance + I(contiguity>0)

fitsur.f.d.l.i <- systemfit(list(coop = coop.f.d.l.i, conf = conf.f.d.l.i), data=final.data.dyad, method = "OLS")
summary(fitsur.f.d.l.i)
## ----------------------------------------------------------------------------


## 3. DY LDV relative events --------------------------------------------------
# look at relative events simultaneously with total environmentIsLocked(

rel.f.d.l.t <- I(log(post.coopTotal + logdelta) - log(abs(post.hostTotal) + logdelta))  ~  log(coopTotal + logdelta) + log(abs(hostTotal) + logdelta) + 
  (KOFEcGIdf.x.0 + KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dep1 * dep2 + dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) +
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2 + 
  alliance + I(contiguity>0)


tot.f.d.l.t <- I(log(post.coopTotal + abs(post.hostTotal) + logdelta))  ~ log(coopTotal + logdelta) + log(abs(hostTotal) + logdelta)  + 
  (KOFEcGIdf.x.0 + KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dep1 * dep2 + dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) +
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2 + 
  alliance + I(contiguity>0)

fitsur.f.d.l.t <- systemfit(list(coop = rel.f.d.l.t, conf = tot.f.d.l.t), data=final.data.dyad, method = "OLS")
summary(fitsur.f.d.l.t)
## ----------------------------------------------------------------------------


## 4.DY Year FE without Trade DEP variables -----------------------------------
coop.f.d.wod <- log(post.coopTotal + logdelta)  ~   (KOFEcGIdf.x.0 + KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) + 
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2 + as.factor(year)+ 
  alliance + I(contiguity>0)

conf.f.d.wod  <- log(abs(post.hostTotal) + logdelta)  ~  (KOFEcGIdf.x.0 + KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) + 
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2 + as.factor(year)+ 
  alliance + I(contiguity>0)

fitsur.f.d.wod  <- systemfit(list(coop = coop.f.d.wod, conf = conf.f.d.wod), data=final.data.dyad, method = "OLS")
summary(fitsur.f.d.wod )
## ----------------------------------------------------------------------------


## 5. DY LDV without Trade DEP variables --------------------------------------
coop.f.d.l.wod  <- log(post.coopTotal + logdelta)  ~  log(coopTotal + logdelta) + log(abs(hostTotal) + logdelta) + (KOFEcGIdf.x.0 + KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) + log(wdi_gdpcappppcon2011.x) + 
  log(wdi_gdpcappppcon2011.y) + dem1 * dem2 + 
  alliance + I(contiguity>0)

conf.f.d.l.wod  <- log(abs(post.hostTotal) + logdelta)  ~ log(coopTotal + logdelta) + log(abs(hostTotal) + logdelta)  + (KOFEcGIdf.x.0 + KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dyadic.balance + caprat +  log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) + log(wdi_gdpcappppcon2011.x) + 
  log(wdi_gdpcappppcon2011.y)+ dem1 * dem2+ 
  alliance + I(contiguity>0) 

fitsur.f.d.l.wod  <- systemfit(list(coop = coop.f.d.l.wod, conf = conf.f.d.l.wod), data=final.data.dyad, method = "OLS")
summary(fitsur.f.d.l.wod )
## ----------------------------------------------------------------------------


## 6. DY LDV with globalization scale between 0 and 1 -------------------------
coop.f.d.l.01 <- log(post.coopTotal + logdelta)  ~  log(coopTotal + logdelta) + log(abs(hostTotal) + logdelta) +
  (KOFEcGIdf.x + KOFEcGIdf.y) * (major.sym + major.asym1 + major.asym2) + 
  dep1 * dep2 + dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) +
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2+ 
  alliance + I(contiguity>0) 

conf.f.d.l.01 <- log(abs(post.hostTotal) + logdelta)  ~ log(coopTotal + logdelta) + log(abs(hostTotal) + logdelta) + 
  (KOFEcGIdf.x + KOFEcGIdf.y) * (major.sym + major.asym1 + major.asym2) + 
  dep1 * dep2 + dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) +
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2 + 
  alliance + I(contiguity>0)

fitsur.f.d.l.01 <- systemfit(list(coop = coop.f.d.l.01, conf = conf.f.d.l.01), data=final.data.dyad, method = "OLS")
summary(fitsur.f.d.l.01)
## ----------------------------------------------------------------------------


## 7. DY LDV with pta  -------------------------------------------------------
coop.f.d.l.pta <- log(post.coopTotal + logdelta)  ~  log(coopTotal + logdelta) + log(abs(hostTotal) + logdelta) +
  (KOFEcGIdf.x.0 + KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dep1 * dep2 + dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) +
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2+ 
  alliance + I(contiguity>0) + dyadpta + log(pta1+1) + log(pta2+1)

conf.f.d.l.pta <- log(abs(post.hostTotal) + logdelta)  ~ log(coopTotal + logdelta) + log(abs(hostTotal) + logdelta) + 
  (KOFEcGIdf.x.0 + KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dep1 * dep2 + dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) +
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2+ 
  alliance + I(contiguity>0) + dyadpta + log(pta1+1) + log(pta2+1)

fitsur.f.d.l.pta <- systemfit(list(coop = coop.f.d.l.pta, conf = conf.f.d.l.pta), data=final.data.dyad, method = "OLS")
summary(fitsur.f.d.l.pta)
## ----------------------------------------------------------------------------

## 8. DY LDV KOF with major count  -------------------------------------------------------
coop.f.d.l.count <- log(post.coopTotal + logdelta)  ~  log(coopTotal + logdelta) + log(abs(hostTotal) + logdelta) +
  (KOFEcGIdf.x.0 + KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dep1 * dep2 + dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) +
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2+ 
  alliance + I(contiguity>0) 

conf.f.d.l.count <- log(abs(post.hostTotal) + logdelta)  ~ log(coopTotal + logdelta) + log(abs(hostTotal) + logdelta) + 
  (KOFEcGIdf.x.0 + KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dep1 * dep2 + dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) +
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2+ 
  alliance + I(contiguity>0) 

final.data.dyad.count <- final.data.dyad %>% 
  mutate(KOFEcGIdf.x.0=log(KOF.count.maj.1+1),
         KOFEcGIdf.y.0=log(KOF.count.maj.2+1))
fitsur.f.d.l.count <- systemfit(list(coop = coop.f.d.l.count, conf = conf.f.d.l.count), data=final.data.dyad.count, method = "OLS")
summary(fitsur.f.d.l.count)
## ----------------------------------------------------------------------------

## 8.b DY LDV KOF with major count (both)  -------------------------------------------------------
coop.f.d.l.countboth <- log(post.coopTotal + logdelta)  ~  log(coopTotal + logdelta) + log(abs(hostTotal) + logdelta) +
  (KOFEcGIdf.x.0 + KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dep1 * dep2 + dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) +
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2+ 
  alliance + I(contiguity>0) 

conf.f.d.l.countboth <- log(abs(post.hostTotal) + logdelta)  ~ log(coopTotal + logdelta) + log(abs(hostTotal) + logdelta) + 
  (KOFEcGIdf.x.0 + KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dep1 * dep2 + dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) +
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2+ 
  alliance + I(contiguity>0) 

final.data.dyad.countboth <- final.data.dyad %>% 
  mutate(KOFEcGIdf.x.0=log(KOF.count.maj.both.1+1),
         KOFEcGIdf.y.0=log(KOF.count.maj.both.2+1))
fitsur.f.d.l.countboth <- systemfit(list(coop = coop.f.d.l.countboth, conf = conf.f.d.l.countboth), data=final.data.dyad.countboth, method = "OLS")
summary(fitsur.f.d.l.countboth)
## ----------------------------------------------------------------------------



## 9. DY LDV with opponent's coop and conflict  -------------------------------------------------------
coop.f.d.l.4iv <- log(post.coopTotal + logdelta)  ~  log(coopTotal + logdelta) + log(abs(hostTotal) + logdelta) +
  log(coopTotal.reverse + logdelta) + log(abs(hostTotal.reverse) + logdelta) +
  (KOFEcGIdf.x.0 + KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dep1 * dep2 + dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) +
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2+ 
  alliance + I(contiguity>0) 

conf.f.d.l.4iv <- log(abs(post.hostTotal) + logdelta)  ~ log(coopTotal + logdelta) + log(abs(hostTotal) + logdelta) + 
  log(coopTotal.reverse + logdelta) + log(abs(hostTotal.reverse) + logdelta) +
  (KOFEcGIdf.x.0 + KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dep1 * dep2 + dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) +
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2+ 
  alliance + I(contiguity>0) 

fitsur.f.d.l.4iv <- systemfit(list(coop = coop.f.d.l.4iv, conf = conf.f.d.l.4iv), data=final.data.dyad, method = "OLS")
summary(fitsur.f.d.l.4iv)
## ----------------------------------------------------------------------------

## 10. DY LDV KOF replaced with dep  -------------------------------------------------------
coop.f.d.l.dep <- log(post.coopTotal + logdelta)  ~  log(coopTotal + logdelta) + log(abs(hostTotal) + logdelta) +
  (KOFEcGIdf.x.0 + KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) +
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2+ 
  alliance + I(contiguity>0) 

conf.f.d.l.dep <- log(abs(post.hostTotal) + logdelta)  ~ log(coopTotal + logdelta) + log(abs(hostTotal) + logdelta) + 
  (KOFEcGIdf.x.0 + KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) +
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2+ 
  alliance + I(contiguity>0) 

final.data.dyad.dep <- final.data.dyad %>% 
  mutate(KOFEcGIdf.x.0=dep1,
         KOFEcGIdf.y.0=dep2)
fitsur.f.d.l.dep <- systemfit(list(coop = coop.f.d.l.dep, conf = conf.f.d.l.dep), data=final.data.dyad.dep, method = "OLS")
summary(fitsur.f.d.l.dep)
## ----------------------------------------------------------------------------

## 11. DY LDV with 3rd party events  -------------------------------------------------------
coop.f.d.l.3rdevents <- log(post.coopTotal + logdelta)  ~  log(coopTotal + logdelta) + log(abs(hostTotal) + logdelta) +
  log(coop.all.x + logdelta) + log(abs(host.all.x) + logdelta) +
  log(coop.all.y + logdelta) + log(abs(host.all.y) + logdelta) +
  (KOFEcGIdf.x.0 + KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dep1 * dep2 + dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) +
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2+ 
  alliance + I(contiguity>0) 

conf.f.d.l.3rdevents <- log(abs(post.hostTotal) + logdelta)  ~ log(coopTotal + logdelta) + log(abs(hostTotal) + logdelta) + 
  log(coop.all.x + logdelta) + log(abs(host.all.x) + logdelta) +
  log(coop.all.y + logdelta) + log(abs(host.all.y) + logdelta) +
  (KOFEcGIdf.x.0 + KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dep1 * dep2 + dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) +
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2+ 
  alliance + I(contiguity>0) 

fitsur.f.d.l.3rdevents <- systemfit(list(coop = coop.f.d.l.3rdevents, conf = conf.f.d.l.3rdevents), data=final.data.dyad, method = "OLS")
summary(fitsur.f.d.l.3rdevents)
## ----------------------------------------------------------------------------

## 12. DY LDV without zero interactions  -------------------------------------------------------
fitsur.f.d.l.wozero <- systemfit(list(coop = coop.f.d.l, conf = conf.f.d.l), 
                          data=final.data.dyad %>% filter(coopTotal!=0|hostTotal!=0),
                          method = "OLS")
summary(fitsur.f.d.l.wozero)
## ----------------------------------------------------------------------------

## 13. DY Mean KOF FE -----------------------------------------------------------------
final.data.dyad <- final.data.dyad %>%
  group_by(year) %>%
  mutate(kofaverage=mean(KOFEcGIdf.x.0, na.rm=TRUE)) %>%
  ungroup()

coop.f.d.kof <- log(post.coopTotal + logdelta)  ~   (KOFEcGIdf.x.0 + KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dep1 * dep2 + dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) + 
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2 + kofaverage+ 
  alliance + I(contiguity>0)

conf.f.d.kof <- log(abs(post.hostTotal) + logdelta)  ~  (KOFEcGIdf.x.0 + KOFEcGIdf.y.0) * (major.sym + major.asym1 + major.asym2) + 
  dep1 * dep2 + dyadic.balance + caprat + log(mindist + 1) + log(wdi_pop.x) + log(wdi_pop.y) + 
  log(wdi_gdppppcon2011.x) + log(wdi_gdppppcon2011.y) + dem1 * dem2 + kofaverage + 
  alliance + I(contiguity>0)

fitsur.f.d.kof <- systemfit(list(coop = coop.f.d.kof, conf = conf.f.d.kof), data=final.data.dyad, method = "OLS")
summary(fitsur.f.d.kof)
